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Abstract 



Regular firing neurons can be seen as oscillators. The phase-response curve (PRC) 
describes how such neurons will respond to small excitatory perturbations. Knowledge 
of the PRC is important as it is associated to the excitability type of neurons and their 
capability to synchronize in networks. In this work we present a novel method to estimate 
the PRC from experimental data. We assume that continuous noise signal can be discretized 
into independent perturbations at evenly spaced phases and predict the next spike based 
on these independent perturbations. The difference between the predicted next spike made 
at every discretized phase and the actual next spike time is used as the error signal used to 
optimize the PRC. We test our method on model data and experimentally obtained data 
and find that the newly developed method is robust and reliable method for the estimation 
of PRCs from experimental data. 

1 Introduction 

The phase response curve (PRC) is a property of oscillators that defines how the phase of the 
next cycle is shifted by an incoming perturbation at a particular phase. In the case of regular 
firing neurons, the PRC describes how the timing of the next spike is shifted. The PRC is 
of interest for several reasons. First, by using the PRC a prediction can be made about the 
synchronizing properties of synoptically coupled neurons. Nonnegative PRCs (i.e., type-I PRCs 
or monophasic PRCs) do not allow for synchronization via excitatory connections while PRCs 
with a negative part (i.e., type-II PRCs or biphasic PRCs) allow for synchronization in networks 
with excitatory coupling and small conduction delays [1] . Second, the PRC type indicates which 
bifurcation type leads to the transition from non-spiking to spiking behavior [5]. Moreover, the 
PRC type correlates with the type of excitability of a neuron [21 II] • As such, the PRC is an 
important characterizing feature of regularly firing neurons. 

The PRC can be easily understood when projecting the higher- dimensional dynamics of a 
model of a neuron to a one-dimensional system, i.e., when a neuron is seen as an oscillator, with 
its state determined by the phase 4>. The unperturbed state of the neuron is then described as 
^ — Lj, with ai = ^ and T the average inter-spike interval. With the introduction of noise 
X (<f>{t)) of small amplitude a, the state of the neuron is described as 



where Z {(p{t)) is the PRC. By using a shorthand notation (f) = (t){t) we can describe the spike 
times of the regularly firing neuron as 
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S = UJ + a ■ X {(j)) ■ Z (0) 
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Thus, the firing times of the neuron is determined by the period w, the noise signal and the PRC 
Z{cl>). 

The PRC of a neuron is usually estimated from experimental data by the direct method 
in which small, brief current-pulses are injected at different phases and the phase shifts of the 
next spike are recorded. By following the definition of the PRC, the phase of the perturbation 
is plotted on the x-axis and the resultant shift in normalized spike time (phase) on the y-axis. 
In practice, neuron spike times display a large amount of jitter and therefore it is necessary to 
measure the spike time shifts in hundreds of inter-spike intervals at random phases to span all 
phases and to obtain a reliable readout. 

An alternative method was proposed by Izhikevich , where continuous current fluctuation 
was injected instead of brief pulses. Subsequently, by summing the effects of the injected fluc- 
tuations between two spikes as predicted by a a candidate PRC, the time of the next spike is 
estimated. The summed difference between the observed and predicted spike times is then used 
as an error signal to optimize the PRC. In practice this method does not converge within a rea- 
sonable number of optimization rounds (e.g., 5000 iterations with a basic simplex optimization 
algorithm). The inability to converge is likely caused by (i) the loss of temporal information 
when the errors are summed, and (ii) the fact that many candidate PRCs can have the same 
error value. We thus developed an extension of this method which does not sum over these two 
dimensions and instead uses an array of error signals that preserve the phase information to 
optimize the PRCs. 

Our novel method, the STEP (Standardized Error Prediction) method, allows for an reliable 
estimation of the PRC with relatively few spikes and at high noise levels. The STEP method 
is a numerical method based on the assumption that the effects of small perturbations on the 
PRCs are independent and treating them independently preserves temporal information. 

2 Methods 

2.1 Experimental procedure 

We consider a neuron that is brought to regular firing by a constant current injection (step 
current Ig)- When stabilized, the considered neuron fires with a period oi uj ~ where T is 
the average interspike interval. A continuous fiuction is imposed on top of Is and the combined 
signal is defined as x{t) (the total current injected into the neuron). In the STEP method, 

n— 3 

the PRC is approximated by a truncated Fourier series z{(j)) « ^ asin(n(/)) 4- bcos(n0) where 

n=0 

z{(j)) is a candidate PRC of which a and b are subject to optimization. The optimization is 
guided by the distance between predicted spike times Si to the real spike times s^. We implicitly 
normalize the duration of the fluctuating signal between subsequent spikes by equating this 
interval to [0, 27r] (s^ = and s^+i = 27r). Then, we discretize the period of the oscillation into 
N bins,{0i, 02, • ■ • , Subsequently, for every phase bin based on the fluctuating signal 

x{t), we compute the predicted next spike time s^j = Si-i -t- x{t{(j)-j)) ■ z{(j)-j) and formulate an 

error signal Ei j = \J [s^ j — Si)"^ . Hence, for every spike we obtain a error signal for each phase 
4>j and use Eij with least-squares optimization to optimize the a and b of the truncated Fourier 
series. Figure [1] illustrates this procedure. (In all shown simulations, we use 100 bins.) 
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Figure 1: Experimental procedure of the STEP method. The input fluctuation between two 
spikes (xn) are scaled to the length of the average ISI. Then, at a predetermined number of 
'bins' in the normalized phase ({^1, (/)2, ■ ■ ■ , 0n}), the next spike time is predicted according to 
a candidate PRC Zc{4>). The difference between the true next spike time and the estimated 
spike time is used as error-signal. Due to the binning of normalized phases, we get a temporally 
structured error signal for each spike at each binned phase which allows for optimization with 
least-squares. 



3 



2.2 Model neuron and experimental data 

We use three data sets to test our method. The first data set contains noise-free model data in 
which a single compartmental model neuron (see below) is perturbed at difi^erent phases. This 
set contains 128 perturbations evenly spaced over [0, 27r]. This data set was used as a benchmark 
to compare the stochastic simulations to. 

The second data set contains modeled data from the same single-compartmental model but 
with an additionally injected fluctuating current. The fluctuations are generated through a sta- 
tionary Orstein-Uhlenbeck process around a given mean value and parametrized by the reversion 
rate {g — 0.1) and 4 different volatility levels {D = le~^, 5e~'*, le~^, 5e~^). The advantage of the 
modeled data (with continuous fluctuations) is that we know the excitability type with certainty 
because small perturbations do not change the PRC type |2] . 

The single compartment model neuron is a modification of the model developed by Golomb 
and Amitai [Sj, as modified in [6]. It uses a Hodgkin Huxley type formalism to model neural 
spiking behavior 

Cm ^ = -m^/igjva (y " ENa) ~ "teDi?, (V - Ek) - S5^s - Er) - gieak ~ Eleak) ~ hn] 

(3) 

and ^ — T {V) {x — x^o (V^))where V is the membrane potential, the maximum conductance 
for ion x andEx the reversal potential for ion x. The parameter values can be found By 
turning the adaptation current on or off, this model switches between type-II or type-I excitability 
[6], respectively. The model is simulated using NEURON [7j and simulates 100s of neuronal time 
which results in approximately 950 spikes. 

The third data set contains experimental data recorded from a layer 2/3 pyramidal cell of 
the mouse cortex with the whole cell patch-clamp techniques in vitro. Standard patch-clamp 
techniques as in [5] were used. Membrane potential voltage data and the injected fluctuations 
were digitized at 40 kHz. Two levels of fluctuation amplitudes (50 pA and lOOpA) were tested 
and resulted in 655 and 647 spikes used; only spikes that satisfied 0.1 x I SI < I SI < 2 x ISI 
were used. 

3 Results 

To calibrate the data, we injected short pulses at different phases and directly recorded the 
phase shift. Afterwards, we ran the STEP method on the perturbation data and found that the 
resulting PRC matches well with the directly determined PRC (Figure [2|. The slight mismatch 
at the end of the period in type-I data and middle of the type-II data is due to the short expansion 
(3rd order) of the Fourier series. Longer expansion tend to provide a better fit but are prone to 
over-fitting (not shown). 

Having established that STEP method works well with simple perturbation, we tested our 
method on the more realistic case in which a continuous fluctuation signal perturbs the neuron. 
The PRCs for the four different noise levels Figure [3] (bottom left) and indicate that at all noise 
levels the PRC correctly displays type-II characteristics (a substantial negative part). 

An important feature of any method to estimate the PRC is its sensitivity to the number of 
spikes because large numbers of spikes are hard to obtain experimentally. For this reason, we 
tested our method with 50, 100 and 500 spikes (Figure [3) right). We observe that the STEP 
algorithm provides a fair estimation with a minimum of 100 spikes. This number of required 
spikes is much lower than most other methods, for instance the STA method [8j is shown to work 
with 7000 spikes, Galan's method [5] with 480 spikes and the WSTA method [TO] with 'a few 
hundred' spikes. 
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Figure 2: Calibration of the STEP method. PRCs obtained with the STEP method are compared 
with those obtained from the same data by using the direct method. 




Figure 3: Performance of the STEP method with in terms of fluctuation level and number of 
spikes. Top panel: the injected fluctuation signals. Bottom left: the PRC generated by our 
method for the four different noise levels. All PRCs are of type-II as they should be. Bottom 
right: the estimated PRC based on 50, 100 and 500 spikes (for the model data at fluctuation level 
1). The PRCs indicate are type-II and from 100 spikes and more the generated PRC converges. 



5 




Figure 4: PRC generated by the STEP method on experimental data. To investigate the reha- 
bihty of the produced PRC we generated two additional PRCs each with only half the available 
spiking data. Since all three PRCs are similar (all 3 per noise level, and all 6 over both noise 
levels) , the resulting PRCs can be considered reliable representation of the spiking data. 

Finally, we tested our method on experimental data obtained from mouse layer 2 /3 pyramidal 
neurons in the visual cortex. Unfortunately, the PRC types for such neurons are not known as 
there is experimental evidence that Layer 2/3 neurons can possess both types of PRC [S]. To give 
an indication that our method works well with experimental data we used a strategy outlined in 
[9] in which we compare the PRC produced with all spike data to two PRCs produced with only 
half of all spiking data. If all the resulting PRCs resemble each other, the PRCs are unlikely 
to contain under- or over-fitting artifacts and can therefore be considered representative for the 
provided data. Figure [4] shows these results. 

We may state that the STEP method is a reliable method to estimate the PRC from exper- 
imental data. The advantages of the STEP method are the low number of required spikes and 
the robustness at higher noise levels. 

The authors thank Drs. Stijn Vanderlooy, Yasuhiro Tsubo and Adam Ponzi for fruitful 
discussions. 
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